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Abstract 

Quantum molecular dynamic (QMD) simulations have been applied to study the thermophysical 
properties of liquid xenon under dynamic compressions. The equation of state (EOS) obtained 
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i-rt . from QMD calculations are corrected according to Saha equation, and contributions from atomic 

Q I ionization, which are of predominance in determining the EOS at high temperature and pressure. 



CJ 



0^ 



are considered. For the pressures below 160 GPa, the necessity in accounting for the atomic 



^ . ionization has been demonstrated by the Hugoniot curve, which shows excellent agreement with 

00 , previous experimental measurements, and three levels of ionization have been proved to be sufficient 



at this stage. 
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I. INTRODUCTION 

Accurate understandings of the thermodynamic properties of materials under extreme 
conditions are of great scientific interest^. The relative high temperatures (several eV) and 
densities (several g/cm^) produce the so-called "warm dense matter" (WDM), which pro- 
vides an active research platform by combining the traditional plasma physics and condensed 
matter physics. WDM usually appears in shock or laser heated solids, inertial confinement 
fusion, and giant planetary interiors, where simultaneous dissociations, ionizations, and de- 
generations make modelling of the dynamical, electrical, and optical properties of WDM 
extremely challenging. 

The study of xenon in the warm dense region has long been focused because of the 
potential applications in power engineering and the significance in astrophysics^-. Various 
experimental measurements^"- and theoretical models^""— have been applied to probe the 
EOS of xenon. Since the ionization equilibrium are not interfered with the dissociation 
equilibrium between molecules and atoms, xenon, which is characterized by monoatomic 
molecule and close shell electronic structure, is particularly suitable for the investigation of 
the high pressure behavior under extreme conditions. Despite its simple atomic structure, 
modelling the EOS of xenon under dynamic compression, especially the Hugoniot curve, is an 
essential issue. For instance, 0-K isotherm, which is important in determining the physical 
properties of solid xenon under static pressure, has been modelled by augmented-plane-wave 
(APW) electron band theory method^. While, the EOS of xenon under dynamic compres- 
sion has been calculated through the fluid perturbation theory (FPT), and compared with 
experimental data^^. On the other hand, the chemical picture of fluid xenon can be simply 
described through the concept of plasma physics, deep insight into the shock compressed 
xenon has been extracted by the self-consistent fluid variational theory (SFVT)— . 

However, the electronic structure has been convinced to be the key role in determining 
the thermophysical properties of WMD^. Due to the intrinsic approximations of these 
classical methods, fully quantum mechanical description of xenon under shock compressions 
still remains to be presented and understood. Meanwhile, QMD simulations, where quantum 
effects are considered by the combinations of classical molecular dynamics for the ions and 
density functional theory (DFT) for electrons, have already been proved to be successful in 
describing WDlsM'^. 



In the present paper, QMD simulations, where free electrons induced by the high tem- 
peratures and pressures are considered, have been used to calculate the thermophysical 
properties of warm dense xenon. The EOS are determined for a wide range of densities 
and temperatures, and compared with experimental measurements and different theoretical 
models. The Hugoniot curves, where different levels of ionization are considered, are then 
derived from the EOS. The rest of the paper is organized as follows. The simulation details 
are briefly described in Sec. Ull Corrections to QMD data and the EOS are presented in Sec. 
mil Finally, we close our paper with a summary of our main results. 

II. COMPUTATIONAL METHOD 

In this work, we have performed simulations for shock compressed xenon by employing 
the Vienna Ab-initio Simulation Package (VASP), which was developed at the Technical Uni- 
versity of Viennar^'^. N atoms in a fixed volume supercell, repeated throughout the space 
periodically, form the elements of simulations. After introducing Born-Oppenheimer ap- 
proximation, the thermo-activated electrons are fully quantum mechanically treated through 
plane-wave, finite-temperature DFT— . The electronic states are populated according to the 
Fermi-Dirac distribution at temperature Tg. The exchange correlation functional is deter- 
mined by generalized gradient approximation (GGA) with the parametrization of Perdew- 
Wang 91^-^. Projector augmented wave (PAW) pseudopotential^^ are used to present the 
ion-electron interactions. The electronic structure are obtained from a series of given ionic 
positions, which are subsequently varied according to the forces calculated within the frame- 
work of DFT via the Hellmann-Feynman theorem at each molecular dynamics step. The 
simulations are executed with the isokinetic ensemble (NVT), where the ionic temperature 
(Tj) is controlled by Nose thermostat^, and the system is kept in local thermo-dynamical 
equilibrium by setting the electron temperature Tg equal to Tj. 

The bulk of the QMD simulations are performed using F point to sample the Brillouin 
zone in molecular dynamic simulations with the plane-wave cutoff of 600.0 eV, because EOS 
(pressure and energy) can only be modified within 5% for the selection of higher number 
of k points and larger cutoff energy. A number of 64 atoms is included in the simulation 
box. To cover typical Hugoniot points with respect to the data obtained from experiments, 
the densities adopted in our simulations range from 2.965 to 10.0 g/cm^ and temperatures 
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FIG. 1: (Color online) Charge density distribution for different conditions: p = 5.0 g/cm^, T = 
1000 K (left panel); and p = 9.5 g/cm^, T = 40000 K (right panel). 

between 165 and 50000 K. Each density and temperature point is typically simulated for 4 
~ 6 ps, and the time step selected for the integrations of atomic motion lies between 0.5 
and 2.0 fs with respect to different conditions. Then, the subsequent 1 ps at equilibrium are 
used to calculate EOS as running averages. 

III. RESULTS AND DISCUSSION 

A. EFFECTS OF FREE ELECTRONS IN WARM DENSE XENON 



Internal energy for shock compressed materials can be expressed as: 

Etot = Eion + Tion + / drVext{r)n{r) + ^ drVH{r)n{r) + E^c + Teiectron + AE/, (1) 

where Eion and Tion correspond to the interaction energy and the kinetic energy of the 
bare nuclear. The third term comes from the fixed external potential Vext{f) (in most cases 
the potential is due to the classical nuclei), in which the electrons move. The fourth term 
is the classical electrostatic energy of the electronic density, which can be obtained from 
the Hartree potential. Exc and Teiectron denote exchange- correlation energy and the kinetic 
energy of electrons. Finally, the ionization energy is added in the last term, which stands for 
the energy to generate free electrons and should be seriously considered in the warm dense 
region. In QMD simulations, the last term is excluded, and the simulated results should be 
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FIG. 2: (Color online) The ionized electron number as a function of the density and temperature, 
corrected. The total internal energy is then expressed as: 

Etot = Eqmd + N 2_^ OiiEi, 

i 

where A^ is the total number of atoms for the present system, and a^ is the i-th level 
ionization degree, Ei is the energy required to remove i electrons from a neutral atom, 
creating an z-level ion. Here, Ei equals to ^^ A e,, and A e, is the i-th ionization energy. 

Free electrons can be introduced at high temperatures and pressures by subsequent ion- 
ization of atoms. The ionization precess can be revealed from the change of the charge 
density distribution with increasing density and temperature, as shown in Fig. [H In lower 
density and temperature region, the localization of electrons represents the atomic fluid. 
While, metallic behaviors of fluid xenon are indicated by delocalization of electrons, which 
can be attributed to the thermo-activation of electronic states. Similar phenomena have 
already been found for another nobel gas helium^i and the expanded fluid alkali metal Li^^. 

Although ionizations of electrons can be qualitative described in the framework of finite- 
temperature DFT, contributions to pressure from noninteracting electrons are still lacking. 
Thus, the total pressure should be corrected as: 
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where m presents the mass of xenon atom. The density and temperature are denoted by p 
and T respectively, and ks stands for Boltzmann constant. Eqmd and Pqmd are obtained 
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FIG. 3: (Color online) Component of plasma xenon as a function of temperature at /> = 6.0 g/cm^. 



from VASP. Since, it is difficult to quantify the ionization degree for i-th level from direct 
QMD results, as an alternative selection, we use Salia equation to verify the degree of each 
ionization level: 

(3) 
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where n, is the density of atoms in the i-th state of ionization, that is with i electrons 
removed, gi is the degeneracy of states for the i-ions, ne is the electron density, A is the 
thermal de Broglie wavelength of an electron, rrif. is the electron mass, and Q is the volume 
of the supercell. Here, ai equals to UiVt/N . Similar approximations have already been used 
to study the EOS of liquid deuterium under shock compressions^^, where soften behavior of 
the Hugoniot curve has been found at high pressure, and the simulated results are consist 
with available experiments. 

In the present work, three levels of ionization have been considered (the relative ionization 
energies are 12.1, 21.2, and 32.1 eV respectively), and the relative electron number (equals 
to ^^ iai) as a function of density and temperature are labelled in Fig. [2l Xenon under 
extreme conditions is considered to be partially ionized plasma, where ionization equilib- 
rium: Xe*^*"^-'"'' ^ Xe*"'"+e, stands for the changes in the electronic structure under high 
temperatures and pressures. Three levels of atomic ionization have been taken into account 
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FIG. 4: (Color online) Hugoniot for shocked liquid xenon are compared with previous 
experimental^i^'^i^ and theoretical result a^^'^^ . Specially, the results, which are corrected according 
to different ionization level, are also provided. 

in the present calculations. As an example, the composition for density of 6.0 g/cm^ is 
plotted in Fig. [31 The fraction of charged ions Xe"*", and Xe^"*" increase to maximum, then 
decrease, and the number of Xe^"*" and electrons increase monotonously with temperature, 
while the number of neutral atom Xe decrease with the increase of temperature. The Hugo- 
niot of shocked compressed liquid xenon can be effectively modified with the consideration 
of atomic ionizations, and the detailed discussion can be found in Sec. IIIIBi 

B. THE EQUATION OF STATE 



Let us turn now to see the Hugoniot by investigating the corrected EOS of warm dense 
xenon. Following Lenosky et al, Beule et al, and Hoist et alM'—, smooth functions are 
used to fit the internal energies and pressures by expansions in terms of density (g/cm^) 
and temperature (10'^ K). The corrected QMD data for internal energy (eV/atom) can be 
expanded as follows: 

2 

E = Y,A{T)p\ (5) 



j=0 



Ai(T) = ajoexp[-( — ■ — -Y] + ais + a^T. 



0-12 



(6) 



The total pressure given in GPa can be similarly expanded as: 

2 

P = Y.B^(T)p^^ (7) 

j=0 

BjiT) = bjoexp[~i^-^y] + 6,3 + b,,T. (8) 

Oj2 

The expansion coefficients Oj^ and 6,^ for E and P (accuracy better than 5%) are summarized 
in Tab. |T]and Tab. [Til respectively. 

TABLE I: Coefficients Uik in expansion for the internal energy E. 



I «iO Ctil 0,i2 «i3 ^44 

57.3976 -17.3142 28.2460 -38.0864 1.7336 

1 11.2512 13.0002 -67.2259 -11.4532 0.0057 

2 -0.7180 41.0407 73.9168 0.6057 0.0060 



TABLE II: Coefficients bn^ in expansion for the total pressure P. 
j bjo bji bj2 bj3 bji 

-38.9090 -15.8020 72.3514 63.7065 -0.0607 

1 0.3299 26.0136 -5.7763 -12.5064 -0.0294 

2 0.7302 1.0722 19.2058 0.9256 0.0441 



A crucial measurement of the EOS data of xenon under shock conditions is the Hugoniot, 
which can be derived from the following equation: 

(Eo - E,) + ^{Vo - V,){Po + Pi) = 0, (9) 

where E is the internal energy, P is the pressure, V is the volume, and the subscripts and 1 
denote the initial and shocked state, respectively. The Hugoniot relation, which describes the 
locus of points in {E, P, l^)-space, follows from conservation of mass, momentum, and energy 
for an isolated system compressed by a pusher at a constant velocity. In our simulations, the 
initial density po is 2.965 g/cm^, the relative internal energy £^0=0.01 eV/atom at To=165 
K. The initial pressure Pq of the start point on the Hugoniot can be approximated to be 
zero compared to high pressure of shocked states. 

8 
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FIG. 5: (Color online) Temperature is labelled as a function of density along the Hugoniot. Previous 
experimental'^'^'^ and theoretica l-^^ i^^ results are also included. 

The Hugoniot curve, where results from our simulations and previous predictions from 
both experiments and theories are labelled for comparison, is shown in Fig. |H In order to 
clarify the influence of atomic ionization, results from direct QMD simulations and those 
corrected by Saha equation are also provided in the present work. Agreements can be found 
between experimental measurements and the uncorrected QMD results only below 70 GPa. 
However, at higher pressures [P > 70 GPa), the uncorrected Hugoniot curve shows an 
abrupt increase, which suggests a stiff behavior, and unphysical results are obtained here 
because of the ignorance of atomic ionization. Considerable softening of the Hugoniot can 
be achieved by the corrected QMD calculations, and reliable results are gradually revealed 
through accounting for the first, second, and third level of ionization. Here, the wide- 
range EOS for shock compressed liquid xenon can be well described by the corrected QMD 
simulations (with the consideration of three levels of ionization) below 160 GPa, where good 
agreements are detected between our simulation results and experiments^'^i^'^. Theoretically, 
no Hugoniot points by previous FPT methocP^ are available to be used to compare with 
experiments at pressures higher than 60 GPa, and the Hugoniot curve from SFVT— seems 
to be inconsistent with the results from Nellis et alA. 

Temperature, which is focused as one of the most important parameters in experiments, is 
difficult to be measured because of the uncertainty in determining the optical-intensity loss 



for ultraviolet part of the spectrum in adiabatic or isentropic shock compressions, especially 
for the temperature exceeding 20000 K. QMD simulations provide powerful tools to predict 
shock temperature. Our calculated Hugoniot temperature is shown in Fig. E] as a function 
of the density. For comparison, the previous simulated (FPT^^ and SFVT^^) and experi- 
mental (by Radousky et al}'^ and Urlin et alA) results are also shown in Fig. O One can 
find that the present calculated temperatures for the direct and corrected QMD simulations 
begin to depart at about 10000 K (the respective density is 6.5 g/cm^), which highlights the 
start point of ionization on the Hugoniot. Discontinuous change in temperature, which is 
similar to the behavior of Hugoniot, has been found in the direct QMD simulation results, 
while, smooth variations of temperature along the Hugoniot are predicted by the corrected 
simulation results. The corrected shock temperatures are accordant with the available ex- 
periments. 

IV. CONCLUSION 

In summary, QMD simulations, where contributions from ionized electrons are taken 
into account, are introduced to investigate the EOS of shock compressed liquid xenon, and 
ionization equilibrium: Xe*^*"^-'"'" ^ Xe^'^+e is described by Saha equation. Here, three levels 
of ionization are considered, and it has been demonstrated that the accuracy of the calculated 
EOS can be effectively improved through the present approximations. The corrected EOS 
data are fitted by smooth functions, from which Hugoniot curve are then derived. Good 
agreements have been achieved between our calculated results and those from experimental 
measurements. 
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